JE prépare ma session de travail
setwd("~/GitHub/Cours_2020_UniGE/Cours_Geneve_6")
Je charge les packages utiles
if(!require("stylo")){
install.packages("stylo")
library(stylo)
}
if(!require("cluster")){
install.packages("cluster")
library(cluster)
}
if(!require("FactoMineR")){
install.packages("FactoMineR")
library(FactoMineR)
}
if(!require("factoextra")){
install.packages("factoextra")
library(factoextra)
}
if(!require("smacof")){
install.packages("smacof")
library(smacof)
}
Loading required package: smacof
Loading required package: plotrix
Registered S3 methods overwritten by 'car':
method from
influence.merMod lme4
cooks.distance.influence.merMod lme4
dfbeta.influence.merMod lme4
dfbetas.influence.merMod lme4
Attaching package: ‘smacof’
The following object is masked from ‘package:base’:
transform
Une classification ascendante hiérarchique.
# Je crée un tableau avec des données au hasard
tableau <- cbind(x=c(1,2,5,8,11),y=c(2,4,1,12,10))
colnames(tableau)<- c("Avoir", "Être")
rownames(tableau) <- c("Racine","Molière","Racine","Scarron","Corneille")
tableau
Avoir Être
Racine 1 2
Molière 2 4
Racine 5 1
Scarron 8 12
Corneille 11 10
Calcul de distance, par défaut eucidienne. on peut utiliser le paramètre method=“manhattan” pour une autre distance
distance_tableau <- dist(tableau)
distance_tableau
Racine Molière Racine Scarron
Molière 2.236068
Racine 4.123106 4.242641
Scarron 12.206556 10.000000 11.401754
Corneille 12.806248 10.816654 10.816654 3.605551
Classification/clustering. Par défaut complete linkage. On peut utiliser le paramètre method=“centroid” ou “ward.D”
hc <- hclust(distance_tableau)
plot(hc)
#On customise le graph
plot(hc, main="Ma première classification", xlab="Euclidian distance, \n Complete-linkage clustering" , sub="" , hang =-1)
On peut faire soi-même un delta de Burrows. Pour cela je z-score mes données
#Je transforme mon tableau en data.frame
tableau_df<-as.data.frame(tableau)
#Je z-score mes deux colonnes avec la fonction scale()
tableau_z_a<-scale(tableau_df$Avoir)
tableau_z_b<-scale(tableau_df$Être)
#Je reforme un tableau, avec les nouvelles données
tableau_z<-cbind(tableau_z_a,tableau_z_b)
tableau_z
[,1] [,2]
[1,] -1.05786348 -0.7724598
[2,] -0.81743996 -0.3659020
[3,] -0.09616941 -0.9757388
[4,] 0.62510115 1.2603292
[5,] 1.34637170 0.8537714
Et je produis un nouveau graph
hc <- hclust(distance_tableau, method="ward.D")
plot(hc, main="Ma première classificatoin", xlab="Burrows Delta, \n Ward's method" , sub="" , hang =-1)
Je charge mon corpus
BOYER_AMOURSJUPITERSEMELE_1666<-paste(scan("corpus/BOYER_AMOURSJUPITERSEMELE_1666.txt", what="character", sep=""),collapse=" ")
Read 17734 items
BOYER_ARISTODEME_1648<-paste(scan("corpus/BOYER_ARISTODEME_1648.txt", what="character", sep=""),collapse=" ")
Read 11106 items
CORNEILLEP_ANDROMEDE_1651<-paste(scan("corpus/CORNEILLEP_ANDROMEDE_1651.txt", what="character", sep=""),collapse=" ")
Read 16121 items
CORNEILLEP_MEDEE_1682<-paste(scan("corpus/CORNEILLEP_MEDEE_1682.txt", what="character", sep=""),collapse=" ")
Read 13705 items
DURYER_DYNAMIS_1653<-paste(scan("corpus/DURYER_DYNAMIS_1653.txt", what="character", sep=""),collapse=" ")
Read 16793 items
DURYER_CLITOPHON<-paste(scan("corpus/DURYER_CLITOPHON_1632.txt", what="character", sep=""),collapse=" ")
Read 15265 items
MOLIERE_AMPHITRYON_1668<-paste(scan("corpus/MOLIERE_AMPHITRYON_1668.txt", what="character", sep=""),collapse=" ")
Read 14152 items
MOLIERE_MISANTHROPE_1667<-paste(scan("corpus/MOLIERE_MISANTHROPE_1667.txt", what="character", sep=""),collapse=" ")
Read 16074 items
RACINE_PHEDRE_1697<-paste(scan("corpus/RACINE_PHEDRE_1697.txt", what="character", sep=""),collapse=" ")
Read 13434 items
RACINE_ESTHER_1689<-paste(scan("corpus/RACINE_ESTHER_1689.txt", what="character", sep=""),collapse=" ")
Read 10460 items
SCARRON_DOMJAPHETDARMENIE_1653<-paste(scan("corpus/SCARRON_DOMJAPHETDARMENIE_1653.txt", what="character", sep=""),collapse=" ")
Read 12945 items
SCARRON_ECOLIERDESALAMANQUE_1655<-paste(scan("corpus/SCARRON_ECOLIERDESALAMANQUE_1655.txt", what="character", sep=""),collapse=" ")
Read 14133 items
ROTROU_HERCULEMOURANT_1636<-paste(scan("corpus/ROTROU_HERCULEMOURANT_1636.txt", what="character", sep=""),collapse=" ")
Read 12561 items
ROTROU_DOMBERNARDDECABRERE_1647<-paste(scan("corpus/ROTROU_DOMBERNARDDECABRERE_1647.txt", what="character", sep=""),collapse=" ")
Read 14546 items
SCUDERY_MORTDECESAR_1637<-paste(scan("corpus/SCUDERY_MORTDECESAR_1637.txt", what="character", sep=""),collapse=" ")
Read 11918 items
SCUDERY_ORANTE_1636<-paste(scan("corpus/SCUDERY_ORANTE_1636.txt", what="character", sep=""),collapse=" ")
Read 12408 items
#Je peux regarder ici un exemple (enlever le # au début de la ligne)
#CORNEILLEP_ANDROMEDE_1651
J’en fais un tableau
#Je crée une liste de mes textes
my.corpus.raw = list(BOYER_AMOURSJUPITERSEMELE_1666,
BOYER_ARISTODEME_1648,
DURYER_DYNAMIS_1653,
DURYER_CLITOPHON,
MOLIERE_AMPHITRYON_1668,
MOLIERE_MISANTHROPE_1667,
SCARRON_DOMJAPHETDARMENIE_1653,
SCARRON_ECOLIERDESALAMANQUE_1655)
#Je tokenise mon texte (espace, apostrophe, retrait de la ponctuation)
my.corpus.clean = lapply(my.corpus.raw, txt.to.words)
#Je compte la fréquence des tokens
complete.word.list = make.frequency.list(my.corpus.clean)
#Je transforme le résultat en table, où les résultats de chaque mot/texte sont alignés
table.of.frequencies=make.table.of.frequencies(my.corpus.clean, complete.word.list, relative = F)
processing 8 text samples
combining frequencies into a table...
#Je donne un nom à mes colonnes
row.names(table.of.frequencies)=c("boyer_agamemnon",
"boyer_amoursJupiter",
"duryer_dynamis",
"duryer_clitophon",
"moliere_amphytrion",
"moliere_misanthrope",
"scarron_domJaphet",
"scarron_ecolierSalamanque")
#Je sauvegarde une copie
write.csv(table.of.frequencies, file = "table.of.frequencies.csv",row.names=TRUE)
#je convertis mes données en dataframe
table.of.frequencies = as.data.frame(read.csv(file="table.of.frequencies.csv", sep = ",", header = TRUE, row.names=1, quote = '\"'))
#J'affiche le résultat
#View(table.of.frequencies)
Ce résultat ne m’arrange pas: je vais avoir besoin d’avoir les occurrences en rang, et les textes en colonne, et donc d’inverser mon tableau
#j'utilise la fonction transpose(?)
table.of.frequencies<-t(table.of.frequencies)
#View(table.of.frequencies)
J’observe la distribution de mon corpus.
summary(table.of.frequencies)
boyer_agamemnon boyer_amoursJupiter duryer_dynamis duryer_clitophon moliere_amphytrion
Min. : 0.000 Min. : 0.000 Min. : 0.000 Min. : 0.00 Min. : 0.000
1st Qu.: 0.000 1st Qu.: 0.000 1st Qu.: 0.000 1st Qu.: 0.00 1st Qu.: 0.000
Median : 0.000 Median : 0.000 Median : 0.000 Median : 0.00 Median : 0.000
Mean : 2.208 Mean : 1.425 Mean : 2.127 Mean : 1.92 Mean : 1.783
3rd Qu.: 1.000 3rd Qu.: 0.000 3rd Qu.: 1.000 3rd Qu.: 1.00 3rd Qu.: 1.000
Max. :599.000 Max. :321.000 Max. :589.000 Max. :537.00 Max. :478.000
moliere_misanthrope scarron_domJaphet scarron_ecolierSalamanque
Min. : 0.000 Min. : 0.000 Min. : 0.000
1st Qu.: 0.000 1st Qu.: 0.000 1st Qu.: 0.000
Median : 0.000 Median : 0.000 Median : 0.000
Mean : 2.018 Mean : 1.686 Mean : 1.828
3rd Qu.: 1.000 3rd Qu.: 1.000 3rd Qu.: 1.000
Max. :567.000 Max. :425.000 Max. :423.000
Je peux synthétiser ce résultat avec une boîte à moustache des fréquences (absolues) des mots du corpus.
100% center
Source: stat4decision
maBoite <- boxplot(colSums(table.of.frequencies), main = "Distribution du nombre de mots par de texte", ylab="Nombre de mots (fréq. absolues)", sub="Corpus: distribution")
Je contrôle l’hétérogénéité de mon corpus
#Je dois regrouper les données par auteur
BOYER = colSums(table.of.frequencies[,grepl('boyer_', colnames(table.of.frequencies))])
CORNEILLE = colSums(table.of.frequencies[,grepl('corneille_', colnames(table.of.frequencies))])
DURYER = colSums(table.of.frequencies[,grepl('duryer_', colnames(table.of.frequencies))])
MOLIERE = colSums(table.of.frequencies[,grepl('moliere_', colnames(table.of.frequencies))])
RACINE = colSums(table.of.frequencies[,grepl('racine_', colnames(table.of.frequencies))])
SCARRON = colSums(table.of.frequencies[,grepl('scarron_', colnames(table.of.frequencies))])
ROTROU = colSums(table.of.frequencies[,grepl('rotrou_', colnames(table.of.frequencies))])
#Et le graphique
boxplot(list(BOYER,MOLIERE, CORNEILLE,DURYER,RACINE,SCARRON, ROTROU), names=c('BOYER','MOLIERE', 'CORNEILLE','DURYER','RACINE','SCARRON','ROTROU'), main="Longueur des textes", sub="Corpus: distribution",ylab="Nombre de mots (fréquences absolues)")
Je transforme cette table des fréquences absolues en une table des fréquences relatives
#Je fais une copie de ma table de fréquences
freqs_rel = table.of.frequencies
#Dans cette copie, pour chacun des mots de chaque colonne, je divise le chiffre trouvé par la somme de la colonne
for(i in 1:ncol(freqs_rel)){
freqs_rel[,i] = freqs_rel[,i]/sum(freqs_rel[,i])
}
head(freqs_rel)
boyer_agamemnon boyer_amoursJupiter duryer_dynamis duryer_clitophon moliere_amphytrion
de 0.03211280 0.02666556 0.02459929 0.03310319 0.03173339
et 0.02728784 0.02475494 0.03278050 0.02404143 0.02701985
je 0.01613681 0.01918923 0.01775378 0.01658242 0.02243909
vous 0.02299898 0.01844160 0.01703028 0.00961657 0.01566753
que 0.01704820 0.01910616 0.02276269 0.02416471 0.02077939
le 0.01715542 0.01702941 0.02075913 0.01713722 0.01427339
moliere_misanthrope scarron_domJaphet scarron_ecolierSalamanque
de 0.03162403 0.02984131 0.02738573
et 0.03150669 0.02260918 0.02732099
je 0.02434874 0.02569864 0.02596141
vous 0.03326684 0.02155596 0.02537874
que 0.02018306 0.01530684 0.02026415
le 0.01513729 0.01987080 0.01573223
Je peux (et même je dois) sélectionner dans cette liste les n plus fréquents
#Je prends les premiers rangs. Changez ce chiffre pour changer les résultats par la suite
freqs_rel_mfw = freqs_rel[1:100,]
Regardons à quoi ressemble notre CAH
CAH = agnes(as.dist(dist.wurzburg(t(freqs_rel_mfw))), method="ward")
CAH_orig_aggloCoeff <- CAH$ac
plot(CAH)
Nettoyons tout cela…
plot(CAH,
main="Cluster analysis",
xlab=paste("Cosine Delta\n Agglomerative coefficient =", CAH_orig_aggloCoeff),
hang = -0.1)
Faisons un peu mieux
#on va colorier les vecteurs en fonction de l'auteur (en pratique, la chaîne de caractère avant l'_underscore_)
labels_color = vector(length = length(CAH$order.lab))
labels_color[grep("boyer", CAH$order.lab)] = "darkblue"
labels_color[grep("duryer", CAH$order.lab)] = "darkgrey"
labels_color[grep("scarron", CAH$order.lab)] = "darkred"
labels_color[grep("moliere", CAH$order.lab)] = "darkgreen"
#je reprends mon graph
plot(CAH)
#Je fais un peu de coloriage
factoextra::fviz_dend(CAH,
ylab=paste("Cosine Delta\n Agglomerative coefficient =", round(CAH_orig_aggloCoeff, digit=3)),
k=2, # je cherche à dégager deux clusters
rect = TRUE,
k_colors = rep("black",2), #je garde la couleur des traits en noir
labels_track_height = 0.7, # Je gère la taille de mon rectangle
label_cols = labels_color, #j'applique mes couleurs
horiz=T) # je mets mon graph à l'horizontal, parce que c'est plus pratique
C’est bien mieux, mais pas encore parfait: le code est est trop long, pas très propre. Cachons tout cela dans un coin. Pour cela, compliquons-nous le travail un peu plus en créant une fonction je que conserve dans un autre fichier appelé functions.R.
Je vais donc charger cette fonction avant de l’utiliser.
source("functions.R")
Quatre paramètres ont été prévus * Le nom des données * La légende à afficher * Le coeffichient d’agglomération (récupéré au passage, un peu plus haut) * une variable pour le paramètre labels_track_height
customPlot(CAH,"Original texts, 100 MFW, Culled @ 0%, Distance: wurzburg",CAH_orig_aggloCoeff,0.7)
Je peux donc reproduire ma CAH simplement avec d’autres configurations. Tentons désormais avec plus de MFW
#on passe de 100 à 1000
freqs_rel_mfw = freqs_rel[1:1000,]
#je recalcule ma CAH
CAH = agnes(as.dist(dist.wurzburg(t(freqs_rel_mfw))), method="ward")
CAH_orig_aggloCoeff <- CAH$ac
#Et désormais, je n'ai plus qu'une simple fonction "maison" pour me produire ma visualisation
customPlot(CAH,"Original texts, 1000 MFW, Culled @ 0%, Distance: wurzburg",CAH_orig_aggloCoeff,0.5)
Tentons avec une autre distance: Ruzicka measure ou Minmax (cf.Koppel, M. and Winter, Y. (2014). Determining if two documents are written by the same author. Journal of the Association for Information Science and Technology, 65: 178–87) dont la robustesse a été démontrée par Kestmont et al. ( Kestemont, M.Stover, J.Koppel, M.Karsdorp, K.Daelemans, W., Authorship Verification with the Ruzicka Metric, DH 2016, Jul 2016, Krakow (Poland)).
100% center
(tf pour term frequency)
#je recalcule ma CAH
CAH = agnes(as.dist(dist.minmax(t(freqs_rel_mfw))), method="ward")
CAH_orig_aggloCoeff <- CAH$ac
#Et désormais, je n'ai plus qu'une simple fonction "maison" pour me produire ma visualisation
customPlot(CAH,"Original texts, 1000 MFW, Culled @ 0%, Distance: Minmax",CAH_orig_aggloCoeff,-0.2)
Et encore une autre, en diminuant de nouveau les MFW
#on repasse de 1000 à 100
freqs_rel_mfw = freqs_rel[1:100,]
#je recalcule ma CAH
CAH = agnes(as.dist(dist.wurzburg(t(freqs_rel_mfw))), method="ward")
CAH_orig_aggloCoeff <- CAH$ac
customPlot(CAH,"Original texts, 100 MFW, Culled @ 0%, Distance: wurzburg",CAH_orig_aggloCoeff,0.7)
Tout cela est très bien. Mais qu’est-ce qui est caché derrière ces clusters? Nous pouvons le savoir. Prenons notre dernière CAH, et regardons le processus de construction des clusters:
CAH_2 = as.hclust(CAH)
plot(CAH_2$height, type="h", ylab="height")
Je remarque que j’ai plusieurs clusters. Tentons de voir ce qu’ils contiennent. Classons par 4 notre corpus:
#Je "coupe en deux" ma CAH, ou plutôt je demande une classification en deux des différents
classes = cutree(CAH, k = "4")
#J'ajoute mes classes aux fréquences que j'ai utilisées
Processed_classes = t(freqs_rel_mfw)
Processed_classes = cbind(as.data.frame(Processed_classes), as.factor(classes))
colnames(Processed_classes)[101] = "Classes"
#J'affiche le résultat
Processed_classes[101]
Classes
boyer_agamemnon 1
boyer_amoursJupiter 1
duryer_dynamis 2
duryer_clitophon 2
moliere_amphytrion 3
moliere_misanthrope 3
scarron_domJaphet 4
scarron_ecolierSalamanque 4
Diminuons à deux le nombre de cluster
#Je "coupe en deux" ma CAH, ou plutôt je demande une classification en deux des différents
classes = cutree(CAH, k = "2")
#J'ajoute mes classes aux fréquences que j'ai utilisées
Processed_classes = t(freqs_rel_mfw)
Processed_classes = cbind(as.data.frame(Processed_classes), as.factor(classes))
colnames(Processed_classes)[101] = "Classes"
#J'affiche le résultat
Processed_classes[101]
Classes
boyer_agamemnon 1
boyer_amoursJupiter 1
duryer_dynamis 1
duryer_clitophon 1
moliere_amphytrion 2
moliere_misanthrope 2
scarron_domJaphet 2
scarron_ecolierSalamanque 2
Je peux aussi regarder les valeurs associées à chacun des mots, et ceux associés à chaque cluster:
myClasses = catdes(Processed_classes, num.var = 101)
myClasses
Link between the cluster variable and the quantitative variables
================================================================
Eta2 P-value
je 0.8855821 0.0004896757
amour 0.8040975 0.0025455982
yeux 0.7601916 0.0047643342
moi 0.7512318 0.0053410502
y 0.7445267 0.0058031317
pas 0.7375033 0.0063162127
être 0.7371581 0.0063422162
n 0.7356576 0.0064561179
si 0.7193183 0.0077904004
l 0.7138624 0.0082757444
ne 0.6959504 0.0100182228
faire 0.6770794 0.0121175393
dieux 0.6744536 0.0124323540
qui 0.6573420 0.0146275157
mes 0.6342385 0.0180111247
en 0.6330403 0.0182005359
veux 0.6181419 0.0206765720
m 0.6095886 0.0222026475
ses 0.5960834 0.0247753466
plus 0.5890029 0.0262070694
suis 0.5741345 0.0294079634
comme 0.5457179 0.0363077270
bien 0.5423877 0.0371876525
me 0.5098090 0.0466527660
point 0.5027205 0.0489296469
Description of each cluster by quantitative variables
=====================================================
$`1`
v.test Mean in category Overall mean sd in category Overall sd p.value
amour 2.372484 0.004818126 0.003380804 0.0008299353 0.0016028752 0.01766891
yeux 2.306803 0.001808802 0.001422762 0.0001256128 0.0004427618 0.02106579
si 2.243931 0.009360696 0.007472486 0.0015781275 0.0022263300 0.02483681
l 2.235405 0.017151679 0.015865132 0.0009675451 0.0015227141 0.02539075
dieux 2.172826 0.003528194 0.001978623 0.0014225437 0.0018868420 0.02979338
qui 2.145086 0.009527600 0.008502894 0.0005062741 0.0012638733 0.03194597
mes 2.107052 0.003316553 0.002528891 0.0007548847 0.0009890388 0.03511304
ses 2.042690 0.002753466 0.001970091 0.0007580942 0.0010146507 0.04108310
plus 2.030522 0.005931663 0.005243258 0.0006308718 0.0008969846 0.04230349
suis -2.004730 0.001103213 0.002013334 0.0005036370 0.0012011368 0.04499194
m -2.065701 0.004543850 0.005795718 0.0008642927 0.0016033934 0.03885672
veux -2.080143 0.001166075 0.001530247 0.0001085254 0.0004631936 0.03751245
en -2.105061 0.010588157 0.012344905 0.0015522742 0.0022079724 0.03528599
faire -2.177052 0.001783577 0.002343016 0.0004299235 0.0006798804 0.02947669
ne -2.207182 0.008047953 0.009833099 0.0010665234 0.0021398561 0.02730134
n -2.269274 0.005291597 0.006619240 0.0008650303 0.0015479011 0.02325169
être -2.271587 0.001702773 0.002472563 0.0003633443 0.0008965862 0.02311147
pas -2.272119 0.004991906 0.006090712 0.0006925025 0.0012794955 0.02307934
y -2.282912 0.001141360 0.001708866 0.0003815097 0.0006577039 0.02243556
moi -2.293169 0.003940317 0.005481745 0.0003737930 0.0017784282 0.02183829
je -2.489794 0.017415563 0.021013767 0.0011821199 0.0038235906 0.01278171
$`2`
v.test Mean in category Overall mean sd in category Overall sd p.value
je 2.489794 0.0246119713 0.021013767 0.0013957565 0.0038235906 0.01278171
moi 2.293169 0.0070231732 0.005481745 0.0011974515 0.0017784282 0.02183829
y 2.282912 0.0022763722 0.001708866 0.0002747232 0.0006577039 0.02243556
pas 2.272119 0.0071895170 0.006090712 0.0006163696 0.0012794955 0.02307934
être 2.271587 0.0032423534 0.002472563 0.0005390367 0.0008965862 0.02311147
n 2.269274 0.0079468818 0.006619240 0.0007200349 0.0015479011 0.02325169
ne 2.207182 0.0116182445 0.009833099 0.0012833566 0.0021398561 0.02730134
faire 2.177052 0.0029024541 0.002343016 0.0003371910 0.0006798804 0.02947669
en 2.105061 0.0141016520 0.012344905 0.0010809285 0.0022079724 0.03528599
veux 2.080143 0.0018944192 0.001530247 0.0003899695 0.0004631936 0.03751245
m 2.065701 0.0070475865 0.005795718 0.0011226719 0.0016033934 0.03885672
suis 2.004730 0.0029234554 0.002013334 0.0009875056 0.0012011368 0.04499194
plus -2.030522 0.0045548539 0.005243258 0.0005131879 0.0008969846 0.04230349
ses -2.042690 0.0011867150 0.001970091 0.0005069224 0.0010146507 0.04108310
mes -2.107052 0.0017412296 0.002528891 0.0003817371 0.0009890388 0.03511304
qui -2.145086 0.0074781878 0.008502894 0.0009156383 0.0012638733 0.03194597
dieux -2.172826 0.0004290517 0.001978623 0.0005425611 0.0018868420 0.02979338
l -2.235405 0.0145785853 0.015865132 0.0006251135 0.0015227141 0.02539075
si -2.243931 0.0055842772 0.007472486 0.0005403116 0.0022263300 0.02483681
yeux -2.306803 0.0010367229 0.001422762 0.0002797223 0.0004427618 0.02106579
amour -2.372484 0.0019434819 0.003380804 0.0005637698 0.0016028752 0.01766891
Je commence par faire une ACP.
ACP = PCA(t(freqs_rel_mfw))
Faisons un peu de ménage, ajoutons de l’information, et tentons de comprendre
#categories
get_categories = vector(length = length(colnames(freqs_rel_mfw)))
get_categories[grep("racine", colnames(freqs_rel_mfw))] = "Racine"
get_categories[grep("moliere", colnames(freqs_rel_mfw))] = "Moliere"
get_categories[grep("duryer", colnames(freqs_rel_mfw))] = "Duryer"
get_categories[grep("boyer", colnames(freqs_rel_mfw))] = "Boyer"
ACP = PCA(t(freqs_rel_mfw))
fviz_pca_ind(ACP, col.ind = get_categories, legend="none")
On voit nettement deux axes, qu doivent bien s’appuyer sur quelque chose. Serait-il possible de faire ressortir les mots qui sont à l’origine de ce partitionnement des données?
Continuons notre nettoyage pour y voir plus clair
#categories
get_categories = vector(length = length(colnames(freqs_rel_mfw)))
get_categories[grep("racine", colnames(freqs_rel_mfw))] = "Racine"
get_categories[grep("moliere", colnames(freqs_rel_mfw))] = "Moliere"
ACP = PCA(t(freqs_rel_mfw))
fviz_pca_ind(ACP, col.ind = get_categories, legend="none")
fviz_pca_var(ACP, col.var="contrib",geom.var = "text", select.var = list(contrib =20))+
scale_color_gradient2(
low="yellow", mid="orange",
high="darkred", midpoint=1.6)+theme_bw()
Nous nous rappelons qu’il existe un pourcentage de réalité retenu par axe, que nous pouvons observer, et donc calculer. Dans ce cas précis, il y en a seulement trois axes, car nous n’avons mis que quatre de textes.
ACP$eig
eigenvalue percentage of variance cumulative percentage of variance
comp 1 30.037372 30.037372 30.03737
comp 2 21.041028 21.041028 51.07840
comp 3 15.055366 15.055366 66.13377
comp 4 11.101088 11.101088 77.23485
comp 5 9.201767 9.201767 86.43662
comp 6 7.811868 7.811868 94.24849
comp 7 5.751511 5.751511 100.00000
Un petit graph montre bien la lente perte de significativité des axes :
barplot(ACP$eig[,1], main="Percentage of variance", names.arg=1:nrow(ACP$eig))
Avec l’algorithme de Barnes-Hut,
library(Rtsne)
maRtsne = Rtsne(t(freqs_rel_mfw), dims = 2, initial_dims = 36, perplexity = 0.5, theta = 0.0, check_duplicates = TRUE, pca = TRUE)
plot(maRtsne$Y)
text(maRtsne$Y[,1], maRtsne$Y[,2], labels = row.names(t(freqs_rel_mfw)), cex=.6)
Attention! Si on relance la même commande, on obtient un résultat différent!
library(Rtsne)
maRtsne = Rtsne(t(freqs_rel_mfw), dims = 2, initial_dims = 36, perplexity = 0.5, theta = 0.0, check_duplicates = TRUE, pca = TRUE)
plot(maRtsne$Y)
text(maRtsne$Y[,1], maRtsne$Y[,2], labels = row.names(t(freqs_rel_mfw)), cex=.6)
#encore une fois
library(Rtsne)
maRtsne = Rtsne(t(freqs_rel_mfw), dims = 2, initial_dims = 36, perplexity = 0.5, theta = 0.0, check_duplicates = TRUE, pca = TRUE)
plot(maRtsne$Y)
text(maRtsne$Y[,1], maRtsne$Y[,2], labels = row.names(t(freqs_rel_mfw)), cex=.6)
#encore une fois
library(Rtsne)
maRtsne = Rtsne(t(freqs_rel_mfw), dims = 2, initial_dims = 36, perplexity = 0.5, theta = 0.0, check_duplicates = TRUE, pca = TRUE)
plot(maRtsne$Y)
text(maRtsne$Y[,1], maRtsne$Y[,2], labels = row.names(t(freqs_rel_mfw)), cex=.6)
On peut faire varier la perpléxité (c’est à dire s’autoriser une distorsion plus grande): on passe à 1.5
maRtsne = Rtsne(t(freqs_rel_mfw), dims = 2, initial_dims = 36, perplexity = 1.5, theta = 0.0, check_duplicates = TRUE, pca = TRUE)
plot(maRtsne$Y)
text(maRtsne$Y[,1], maRtsne$Y[,2], labels = row.names(t(freqs_rel_mfw)), cex=.6)
#Encore une fois
maRtsne = Rtsne(t(freqs_rel_mfw), dims = 2, initial_dims = 36, perplexity = 1.5, theta = 0.0, check_duplicates = TRUE, pca = TRUE)
plot(maRtsne$Y)
text(maRtsne$Y[,1], maRtsne$Y[,2], labels = row.names(t(freqs_rel_mfw)), cex=.6)
#Encore une fois
maRtsne = Rtsne(t(freqs_rel_mfw), dims = 2, initial_dims = 36, perplexity = 1.5, theta = 0.0, check_duplicates = TRUE, pca = TRUE)
plot(maRtsne$Y)
text(maRtsne$Y[,1], maRtsne$Y[,2], labels = row.names(t(freqs_rel_mfw)), cex=.6)
#Encore une fois
maRtsne = Rtsne(t(freqs_rel_mfw), dims = 2, initial_dims = 36, perplexity = 1.5, theta = 0.0, check_duplicates = TRUE, pca = TRUE)
plot(maRtsne$Y)
text(maRtsne$Y[,1], maRtsne$Y[,2], labels = row.names(t(freqs_rel_mfw)), cex=.6)
On peut faire varier la perpléxité: on passe à 2
maRtsne = Rtsne(t(freqs_rel_mfw), dims = 2, initial_dims = 36, perplexity = 2, theta = 0.0, check_duplicates = TRUE, pca = TRUE)
plot(maRtsne$Y)
text(maRtsne$Y[,1], maRtsne$Y[,2], labels = row.names(t(freqs_rel_mfw)), cex=.6)
#Encore une fois
maRtsne = Rtsne(t(freqs_rel_mfw), dims = 2, initial_dims = 36, perplexity = 2, theta = 0.0, check_duplicates = TRUE, pca = TRUE)
plot(maRtsne$Y)
text(maRtsne$Y[,1], maRtsne$Y[,2], labels = row.names(t(freqs_rel_mfw)), cex=.6)
#Encore une fois
maRtsne = Rtsne(t(freqs_rel_mfw), dims = 2, initial_dims = 36, perplexity = 2, theta = 0.0, check_duplicates = TRUE, pca = TRUE)
plot(maRtsne$Y)
text(maRtsne$Y[,1], maRtsne$Y[,2], labels = row.names(t(freqs_rel_mfw)), cex=.6)
#Encore une fois
maRtsne = Rtsne(t(freqs_rel_mfw), dims = 2, initial_dims = 36, perplexity = 2, theta = 0.0, check_duplicates = TRUE, pca = TRUE)
plot(maRtsne$Y)
text(maRtsne$Y[,1], maRtsne$Y[,2], labels = row.names(t(freqs_rel_mfw)), cex=.6)
Le Multidimensional scaling est une compression des données en deux dimensions. Il en existe trois types: * classique * métrique * non métrique
Le MDS classique se base sur le calcul de distance, dont il va essayer de préserver l’essentiel.
#Je crée mes données
fit = cmdscale(dist(t(freqs_rel_mfw), method = "manhattan"), eig=TRUE, k=2) # k est le nombre de dimensions souhaité
#Je dessine mon résultat
x = fit$points[,1]
y = fit$points[,2]
plot(x, y, xlab=paste("Coordinate 1 (GOF: ", round(fit$GOF[1] * 100, digits=2), "%)"), ylab=paste("Coordinate 2 (GOF: ", round(fit$GOF[2] * 100, digits=2), "%)"), main="PMD métrique")
text(x, y, labels = row.names(t(freqs_rel_mfw)), cex=.7)
Le MDS métrique, dit aussi ordinal, ne s’intéresse pas à la mesure de distance, mais sa valeur en relation avec les autres paires
Rappelons que le stress permet de mesurer la distortion introduite dans le résultat.
#J'affiche le stress avec le paramètre verbose=T(RUE)
MDSmetrique = mds(dist(t(freqs_rel_mfw), method = "manhattan"), ndim=2, type="ordinal",verbose=T)
Iteration: 1 Stress (raw): 0.01510422 Difference: 0.06233185
Iteration: 2 Stress (raw): 0.00887573 Difference: 0.00622848
Iteration: 3 Stress (raw): 0.00656184 Difference: 0.00231389
Iteration: 4 Stress (raw): 0.00541200 Difference: 0.00114984
Iteration: 5 Stress (raw): 0.00475869 Difference: 0.00065331
Iteration: 6 Stress (raw): 0.00439233 Difference: 0.00036636
Iteration: 7 Stress (raw): 0.00417426 Difference: 0.00021807
Iteration: 8 Stress (raw): 0.00403217 Difference: 0.00014209
Iteration: 9 Stress (raw): 0.00393218 Difference: 0.00009999
Iteration: 10 Stress (raw): 0.00385734 Difference: 0.00007485
Iteration: 11 Stress (raw): 0.00379958 Difference: 0.00005776
Iteration: 12 Stress (raw): 0.00375696 Difference: 0.00004262
Iteration: 13 Stress (raw): 0.00372424 Difference: 0.00003271
Iteration: 14 Stress (raw): 0.00369790 Difference: 0.00002634
Iteration: 15 Stress (raw): 0.00367617 Difference: 0.00002173
Iteration: 16 Stress (raw): 0.00365798 Difference: 0.00001819
Iteration: 17 Stress (raw): 0.00364262 Difference: 0.00001537
Iteration: 18 Stress (raw): 0.00362954 Difference: 0.00001308
Iteration: 19 Stress (raw): 0.00361836 Difference: 0.00001119
Iteration: 20 Stress (raw): 0.00360875 Difference: 0.00000961
Iteration: 21 Stress (raw): 0.00360047 Difference: 0.00000828
Iteration: 22 Stress (raw): 0.00359331 Difference: 0.00000716
Iteration: 23 Stress (raw): 0.00358711 Difference: 0.00000620
Iteration: 24 Stress (raw): 0.00358173 Difference: 0.00000538
Iteration: 25 Stress (raw): 0.00357706 Difference: 0.00000467
Iteration: 26 Stress (raw): 0.00357299 Difference: 0.00000407
Iteration: 27 Stress (raw): 0.00356944 Difference: 0.00000355
Iteration: 28 Stress (raw): 0.00356633 Difference: 0.00000311
Iteration: 29 Stress (raw): 0.00356361 Difference: 0.00000272
Iteration: 30 Stress (raw): 0.00356122 Difference: 0.00000238
Iteration: 31 Stress (raw): 0.00355913 Difference: 0.00000209
Iteration: 32 Stress (raw): 0.00355730 Difference: 0.00000184
Iteration: 33 Stress (raw): 0.00355568 Difference: 0.00000161
Iteration: 34 Stress (raw): 0.00355427 Difference: 0.00000142
Iteration: 35 Stress (raw): 0.00355302 Difference: 0.00000125
Iteration: 36 Stress (raw): 0.00355193 Difference: 0.00000110
Iteration: 37 Stress (raw): 0.00355096 Difference: 0.00000096
plot(MDSmetrique, sub=paste("Stress, ", round(MDSmetrique$stress, digits=2)))
#J'affiche le stress avec le paramètre verbose=T(RUE)
MDSnonMetrique = mds(dist(t(freqs_rel_mfw), method = "euclid"), ndim=2, type="interval",verbose=T)
Iteration: 1 Stress (raw): 0.02349141 Difference: 0.05037561
Iteration: 2 Stress (raw): 0.01787966 Difference: 0.00561176
Iteration: 3 Stress (raw): 0.01617650 Difference: 0.00170315
Iteration: 4 Stress (raw): 0.01541740 Difference: 0.00075910
Iteration: 5 Stress (raw): 0.01499855 Difference: 0.00041885
Iteration: 6 Stress (raw): 0.01474414 Difference: 0.00025441
Iteration: 7 Stress (raw): 0.01458330 Difference: 0.00016085
Iteration: 8 Stress (raw): 0.01447937 Difference: 0.00010393
Iteration: 9 Stress (raw): 0.01441109 Difference: 0.00006828
Iteration: 10 Stress (raw): 0.01436553 Difference: 0.00004556
Iteration: 11 Stress (raw): 0.01433468 Difference: 0.00003085
Iteration: 12 Stress (raw): 0.01431351 Difference: 0.00002117
Iteration: 13 Stress (raw): 0.01429880 Difference: 0.00001471
Iteration: 14 Stress (raw): 0.01428848 Difference: 0.00001033
Iteration: 15 Stress (raw): 0.01428115 Difference: 0.00000732
Iteration: 16 Stress (raw): 0.01427593 Difference: 0.00000523
Iteration: 17 Stress (raw): 0.01427217 Difference: 0.00000376
Iteration: 18 Stress (raw): 0.01426945 Difference: 0.00000272
Iteration: 19 Stress (raw): 0.01426747 Difference: 0.00000197
Iteration: 20 Stress (raw): 0.01426603 Difference: 0.00000144
Iteration: 21 Stress (raw): 0.01426498 Difference: 0.00000105
Iteration: 22 Stress (raw): 0.01426421 Difference: 0.00000077
plot(MDSnonMetrique, sub=paste("Stress, ", round(MDSmetrique$stress, digits=2)))